PS216-PS218_roots-shoots
zagor
2025-10-17
## used (Mb) gc trigger (Mb) max used (Mb)
## Ncells 555854 29.7 1244925 66.5 686457 36.7
## Vcells 1025194 7.9 8388608 64.0 1877160 14.42 palette
my_dir = file.path("..", "reports", "Table7")
if (!dir.exists(my_dir)) {
dir.create(my_dir)
}
fr = file.path('..', 'output', 'Table7_stat.txt')
file.create(fr)## [1] TRUE3 f-ctions
3.1 summary stat
summary_stats <- function(df, measurevar, groupvars, conf.level = 0.95) {
df %>%
dplyr::group_by(across(all_of(groupvars))) %>%
dplyr::summarise(
N = sum(!is.na(.data[[measurevar]])),
mean = mean(.data[[measurevar]], na.rm = TRUE),
sd = sd(.data[[measurevar]], na.rm = TRUE),
.groups = "drop"
) %>%
dplyr::mutate(
se = sd / sqrt(N),
ci = se * qt(conf.level/2 + 0.5, N - 1)
)
}3.2 scale to the average of the reference group
process_data <- function(df, groupvars, measurevar, scale_treatment = "noninoculated") {
df.long = df %>%
tidyr::pivot_longer(cols = 3:ncol(df), names_to = "transcript",
values_to = measurevar, values_drop_na = TRUE)
data.SE = summary_stats(df.long, measurevar, groupvars)
scale_reference = data.SE %>%
dplyr::filter(Treatment == scale_treatment) %>%
dplyr::select(Tissue, transcript, mean) %>%
dplyr::rename(scale_mean = mean)
df.scaled <- dplyr::left_join(df.long, scale_reference, by = c("Tissue", "transcript")) %>%
dplyr::mutate(scaled = .data[[measurevar]] / scale_mean) %>%
dplyr::select(-scale_mean)
df.scaled = df.scaled %>%
dplyr::arrange(dplyr::desc(transcript), dplyr::desc(Treatment))
shoots = df.scaled %>% dplyr::filter(Tissue == "shoots")
roots = df.scaled %>% dplyr::filter(Tissue == "roots")
return(list(
shoots = shoots,
roots = roots
))
}3.3 permutation t-test
- for non-i.i.d, hypotheses
- a non-parametric robust alternative to traditional t-tests
- small sample sizes
- robust to unequal variances
perm_test_by_transcript <- function(mydata.long, measurevar = "measurement", groupvar = "Treatment") {
set.seed(123456)
temp = data.frame()
transcript_levels = levels(mydata.long$transcript)
treatment_levels = levels(mydata.long[[groupvar]])
treatment_pairs = combn(treatment_levels, 2, simplify = FALSE)
for(i in transcript_levels) {
data = mydata.long[mydata.long$transcript == i, ]
k = 8 # nrow(data)
pvalues = purrr::map_dbl(treatment_pairs, function(x) {
subset_data = base::subset(data, data[[groupvar]] %in% x)
res = MKinfer::perm.t.test(
formula = stats::as.formula(base::paste(measurevar, "~", groupvar)),
data = subset_data,
alternative = "two.sided",
mu = 0,
paired = FALSE,
var.equal = FALSE,
conf.level = 0.95,
perm.conf.int = 0.95,
symmetric = TRUE,
p.adjust.method = "BH",
detailed = TRUE,
# split k observations into two nontrivial groups
R = sum(choose(k, 1:(k-1))), # alternative (k, k/2)
set.seed = 123456
)
res$perm.p.value
})
tmp = base::as.data.frame(base::t(pvalues))
colnames(tmp) = purrr::map_chr(treatment_pairs, ~base::paste(.x, collapse = ' vs '))
rownames(tmp) = i
temp <- base::rbind(temp, tmp)
}
return(temp)
}3.4 plot
plot_gene_expression <- function(data.SE
, data.long
, stat.test.sig
, transcripts_excl = c('13-LOX', 'PTI5')
, color_values
, facet_cols = 6
, y_scales = NULL
, dodge_width = 0.8,
plot_title = ""
) {
dodge = position_dodge(width = dodge_width)
data.SE.filtered = dplyr::filter(data.SE, !transcript %in% transcripts_excl)
data.long.filtered = dplyr::filter(data.long, !transcript %in% transcripts_excl)
stat.test.filtered = NULL
if (nrow(stat.test.sig) > 0) stat.test.filtered = dplyr::filter(stat.test.sig, !transcript %in% transcripts_excl)
p = ggplot(data.SE.filtered, aes(x = Treatment, y = mean)) +
geom_point(size=3.5, shape = 22, position = dodge, aes(fill = Treatment), colour = "black") +
geom_point(data = data.long.filtered, aes(x = Treatment, y = scaled, fill = Treatment), colour = "black",
size = 2.0, shape = 21,
position = dodge) +
geom_errorbar(aes(ymin = mean - se, ymax = mean + se), width = 0.3, lwd = 0.5,
position = dodge) +
facet_wrap(~ transcript, ncol = facet_cols, scales = "free", drop = TRUE)
if (!is.null(y_scales)) {
p = p + ggh4x::facetted_pos_scales(y = y_scales)
}
p = p +
# geom_hline(yintercept = 2.5, alpha = 0.0) +
geom_hline(yintercept = 0, alpha = 0.0) +
geom_hline(yintercept = 1, alpha = 0.5, linetype = "dotted", col = "gray45") +
ggtitle(plot_title) +
theme_bw() +
scale_colour_manual(name = "", values = rev(color_values)) +
scale_fill_manual(name = "", values = rev(color_values)) +
labs(x = "", y = "Relative gene expression (+/- SE)"
# , subtitle = "Permutation t-test measurements"
) +
theme(plot.subtitle = element_text(size = 10),
axis.text = element_text(size = 12.5),
axis.text.x = element_text(size = 12.5, angle = 90),
axis.title = element_text(size = 12.5, face = "bold"),
strip.text = element_text(size = 12.5),
title = element_text(size = 15, face = "bold"),
# axis.ticks.x = element_blank(),
legend.key.height = unit(1.5, "cm"),
legend.key.width = unit(1.75, "cm"),
legend.text = element_text(size = 12.5),
legend.title = element_text(size = 12.5),
legend.background = element_rect(fill = "transparent", size = 0.5, linetype = "dotted"),
legend.position = "top",
legend.justification = "right",
plot.background = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
axis.title.y.right = element_blank(),
axis.text.y.right = element_blank(),
# axis.ticks.y = element_blank(),
axis.text.y = element_text(size = 12.5, margin = margin(r = 0)),
panel.spacing = unit(1, "lines"),
strip.background = element_rect(size = 0.5, fill = "transparent", color = NA) ,
panel.border = element_blank(),
axis.line = element_line(color = "black")
)
p = p + theme(legend.position = "none")
if (!is.null(stat.test.filtered) && nrow(stat.test.filtered) > 0) {
p = p + ggpubr::stat_pvalue_manual(
stat.test.filtered,
label = "p.adj.signif",
xmin = "xmin",
xmax = "xmax",
y.position = "y.position",
hide.ns = FALSE,
tip.length = 0.01,
step.increase = 0.0,
inherit.aes = FALSE
)
}
return(p)
}3.5 distribution, variance and effects
- on qPCR data
- quantGenius - quantification of qPCR data using standard curve
- MIQE 2.0: Revision of the Minimum Information for Publication of Quantitative Real-Time PCR Experiments Guidelines
- quantification based on a standard curve inherently involves a log transformation of the input data
- Denote: LOQ values skew distributions
- Denote: distribution tests: low power for small sample size, see QQ of residuals
- if most or all points fall inside the shaded confidence band, the sample’s distribution does not show strong evidence of departure from normality at the plotted sample size and confidence level
- a few isolated points outside the band at the extremes are common with small samples and do not necessarily indicate a severe problem
- Denote: heteroscedasticity: some tests are meant to be used with normally distributed data, but can tolerate relatively low deviation from normality; Levene’s test with mean, more robust test Brown-Forsythe with median, Fligner-Killeen when data are non-normally distributed or when problems related to outliers in the dataset cannot be resolved
- Denote: effect size: Cohen’s d - parametric with assumptions, Wilcoxon Effect Size - the non-parametric alternative; different sensitivity to outliers
dens_and_effect <- function(mydata.long, p_colors) {
# Plot density with x="measurement"
p1 = ggpubr::ggdensity(mydata.long, x = "measurement",
add = "mean", rug = TRUE,
color = "Treatment", fill = "Treatment",
facet.by = 'transcript') +
scale_fill_manual(name = "", values = rev(p_colors)) +
scale_color_manual(name = "", values = rev(p_colors)) +
ggtitle("measurement") +
facet_wrap(~transcript, ncol = 4, scales = "free")
print(p1)
# optional: Plot density with x="scaled"- sam e shape wih a shift
# p2 = ggpubr::ggdensity(mydata.long, x = "scaled",
# add = "mean", rug = TRUE,
# color = "Treatment", fill = "Treatment",
# facet.by = 'transcript') +
# scale_fill_manual(name = "", values = rev(p_colors)) +
# scale_color_manual(name = "", values = rev(p_colors)) +
# ggtitle("scaled") +
# facet_wrap(~transcript, ncol = 4, scales = "free")
# print(p2)
cat(crayon::red('#### #### \nDistribution tests\n#### #### \n'))
distr_tests = mydata.long %>%
dplyr::group_by(transcript) %>%
dplyr::group_map(~ {
x = na.omit(.x$measurement)
if (length(x) < 3) {
return(tibble::tibble(
Shapiro_Wilk = NA_real_,
Anderson_Darling = NA_real_,
Lilliefors_KS = NA_real_,
Jarque_Bera = NA_real_,
DAgostino_Skewness = NA_real_,
n = length(x),
transcript = .y$transcript
))
}
tibble::tibble(
Shapiro_Wilk = tryCatch(shapiro.test(x)$p.value, error = function(e) NA_real_),
Anderson_Darling = tryCatch(nortest::ad.test(x)$p.value, error = function(e) NA_real_),
Lilliefors_KS = tryCatch(nortest::lillie.test(x)$p.value, error = function(e) NA_real_),
Jarque_Bera = tryCatch(tseries::jarque.bera.test(x)$p.value, error = function(e) NA_real_),
DAgostino_Skewness = tryCatch(moments::agostino.test(x)$p.value, error = function(e) NA_real_),
n = length(x),
transcript = .y$transcript
)
}) %>%
dplyr::bind_rows() %>%
dplyr::ungroup() %>%
dplyr::mutate(
Shapiro_Wilk_BH = stats::p.adjust(Shapiro_Wilk, method = "BH"),
Anderson_Darling_BH = stats::p.adjust(Anderson_Darling, method = "BH"),
Lilliefors_KS_BH = stats::p.adjust(Lilliefors_KS, method = "BH"),
Jarque_Bera_BH = stats::p.adjust(Jarque_Bera, method = "BH"),
DAgostino_Skewness_BH = stats::p.adjust(DAgostino_Skewness, method = "BH")
)
combined_results = distr_tests[grep("^n$|transcript|_BH", colnames(distr_tests))]
print(combined_results)
cat(crayon::red("#### #### \nQuantile-Quantile plots\n#### #### \n"))
qq_plot_list = mydata.long %>%
dplyr::group_by(transcript) %>%
dplyr::group_map(~ {
ggpubr::ggqqplot(.x$measurement) +
ggplot2::ggtitle(paste("Raw:", .y$transcript))
})
resid_qq_list = mydata.long %>%
dplyr::group_by(transcript) %>%
dplyr::group_map(~ {
model = lm(measurement ~ Treatment, data = .x)
resids = resid(model)
ggpubr::ggqqplot(resids) +
ggplot2::ggtitle(paste("Residuals:", .y$transcript))
})
resid_density_list = mydata.long %>%
dplyr::group_by(transcript) %>%
dplyr::group_map(~ {
model = lm(measurement ~ Treatment, data = .x)
resids = scale(resid(model)) # standardize residuals
df = data.frame(resids = as.numeric(resids))
ggplot(df, aes(x = resids)) +
geom_density(fill = "steelblue", alpha = 0.6) +
stat_function(fun = dnorm, color = "red", linetype = "dashed") +
ggtitle(paste("Residuals:", .y$transcript)) +
theme_minimal()
})
# Arrange and print all plots
qq_arranged = ggpubr::ggarrange(plotlist = qq_plot_list, ncol = 4, nrow = ceiling(length(qq_plot_list) / 4))
resid_qq_arranged = ggpubr::ggarrange(plotlist = resid_qq_list, ncol = 4, nrow = ceiling(length(resid_qq_list) / 4))
resid_hist_arranged = ggpubr::ggarrange(plotlist = resid_density_list, ncol = 4, nrow = ceiling(length(resid_density_list) / 4))
print(qq_arranged)
print(resid_qq_arranged)
print(resid_hist_arranged)
# Return invisibly
invisible(list(qq = qq_arranged, resid_qq = resid_qq_arranged, resid_hist = resid_hist_arranged))
cat(crayon::red("#### #### \nTest for homogeneity of variance across groups\n#### #### \n"))
cat(crayon::red("#### #### \nLevene\n#### #### \n"))
# Levene's test on measurement
lev_meas = mydata.long %>%
dplyr::group_by(transcript) %>%
rstatix::levene_test(measurement ~ Treatment)
print(lev_meas)
cat(crayon::red("#### #### \nBrown-Forsythe\n#### #### \n"))
# center = "median" switches Levene’s test to the Brown-Forsythe variant
bf_meas = mydata.long %>%
dplyr::group_by(transcript) %>%
rstatix::levene_test(measurement ~ Treatment, center = "median")
print(bf_meas)
cat(crayon::red("#### #### \nFligner\n#### #### \n"))
fk_meas = mydata.long %>%
dplyr::group_by(transcript) %>%
dplyr::group_modify(~ broom::tidy(stats::fligner.test(measurement ~ Treatment, data = .x)))
print(fk_meas)
cat(crayon::red("#### #### \nWilcoxon effect size\n#### #### \n"))
# Wilcoxon effect size on measurement
eff_meas = mydata.long %>%
dplyr::group_by(transcript) %>%
rstatix::wilcox_effsize(measurement ~ Treatment)
print(eff_meas)
# Wilcoxon effect size on scaled
# eff_scaled = mydata.long %>%
# dplyr::group_by(transcript) %>%
# rstatix::wilcox_effsize(scaled ~ Treatment)
# print(eff_scaled)
cat(crayon::red("#### #### \nCohen's d Measure of Effect Size\n#### #### \n"))
# Cohen's d on measurement
coh_meas = mydata.long %>%
dplyr::group_by(transcript) %>%
rstatix::cohens_d(measurement ~ Treatment, paired = FALSE)
print(coh_meas)
invisible(list(
p_measurement = p1,
distr_measurement = combined_results,
levene_measurement = lev_meas,
brownforsythe_measurement = bf_meas,
fligner_measurement = fk_meas,
wilcoxon_measurement = eff_meas,
cohensd_measurement = coh_meas
))
}3.6 ggplot limits
make_scale_params <- function(maxval) {
upper = ifelse(maxval > 2, maxval + 2, maxval + 0.1) # 0.5
step = dplyr::case_when(
maxval >= 50 ~ 25,
maxval >= 30 ~ 10,
maxval >= 15 ~ 5,
maxval >= 10 ~ 2,
maxval > 2 ~ 1,
TRUE ~ 0.5
)
upper_round = ceiling(upper / step) * step
breaks = seq(0, upper_round, by = step)
list(limits = c(0, upper_round), breaks = breaks)
}
build_y_scales_for <- function(names_vec, max_per_transcript) {
max_per_transcript %>%
dplyr::filter(transcript %in% names_vec) %>%
purrr::pmap(function(transcript, max_scaled) {
params = make_scale_params(max_scaled)
lhs = rlang::expr(transcript == !!transcript)
rhs = rlang::expr(ggplot2::scale_y_continuous(limits = !!params$limits, breaks = !!params$breaks))
rlang::new_formula(lhs, rhs)
})
}3.7 test and plot
test_and_plot <- function(data_long_raw
, myorder
, pal
, what
, plot_gene_expression_func
, groupvars = c("Treatment", "transcript")
, y_scales1
, y_scales2
) {
mydata.long = dplyr::as_tibble(data.table::data.table(data_long_raw))
mydata.long$transcript = factor(mydata.long$transcript, levels = myorder)
mydata.long = mydata.long %>% dplyr::arrange(factor(transcript, levels = myorder))
mydata.long.SE = summary_stats(mydata.long, measurevar = "scaled", groupvars = groupvars)
results = dens_and_effect(mydata.long, p_colors = pal)
cat(what, file = fr, append = TRUE, sep = "\n")
cat("\nDistribution tests", file = fr, append = TRUE, sep = "\n")
header = paste(colnames(results$distr_measurement), collapse = "\t")
cat(header, file = fr, append = TRUE, sep = "\n")
output_text = apply(results$distr_measurement, 1, function(row) paste(row, collapse = "\t"))
cat(output_text, file = fr, append = TRUE, sep = "\n")
cat("\nLevene's test for homogeneity of variance across groups", file = fr, append = TRUE, sep = "\n")
header = paste(colnames(results$levene_measurement), collapse = "\t")
cat(header, file = fr, append = TRUE, sep = "\n")
output_text = apply(results$levene_measurement, 1, function(row) paste(row, collapse = "\t"))
cat(output_text, file = fr, append = TRUE, sep = "\n")
cat("\nBrown-Forsythe (Levene with median)", file = fr, append = TRUE, sep = "\n")
header = paste(colnames(results$brownforsythe_measurement), collapse = "\t")
cat(header, file = fr, append = TRUE, sep = "\n")
output_text = apply(results$brownforsythe_measurement, 1, function(row) paste(row, collapse = "\t"))
cat(output_text, file = fr, append = TRUE, sep = "\n")
cat("\nFligner-Killeen test of homogeneity of variances", file = fr, append = TRUE, sep = "\n")
header = paste(colnames(results$fligner_measurement), collapse = "\t")
cat(header, file = fr, append = TRUE, sep = "\n")
output_text = apply(results$fligner_measurement, 1, function(row) paste(row, collapse = "\t"))
cat(output_text, file = fr, append = TRUE, sep = "\n")
cat("\nWilcoxon effect size", file = fr, append = TRUE, sep = "\n")
header = paste(colnames(results$wilcoxon_measurement), collapse = "\t")
cat(header, file = fr, append = TRUE, sep = "\n")
output_text = apply(results$wilcoxon_measurement, 1, function(row) paste(row, collapse = "\t"))
cat(output_text, file = fr, append = TRUE, sep = "\n")
cat("\nCohen's d Measure of Effect Size", file = fr, append = TRUE, sep = "\n")
header = paste(colnames(results$cohensd_measurement), collapse = "\t")
cat(header, file = fr, append = TRUE, sep = "\n")
output_text = apply(results$cohensd_measurement, 1, function(row) paste(row, collapse = "\t"))
cat(output_text, file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
temp = perm_test_by_transcript(mydata.long, measurevar = "measurement", groupvar = "Treatment")
temp$transcript = rownames(temp)
perm = tidyr::gather(temp, contrast, perm.p, colnames(temp)[1]:colnames(temp)[ncol(temp)-1], factor_key=TRUE)
perm$group1 = gsub(' vs.*', '', perm$contrast)
perm$group2 = gsub('.* vs ', '', perm$contrast)
perm$perm.p.adj = p.adjust(perm$perm.p, method = 'BH')
perm$perm.p.adj.signif = 'ns'
perm$perm.p.adj.signif[perm$perm.p.adj < 0.0001] = '**'
perm$perm.p.adj.signif[perm$perm.p.adj < 0.001] = '**'
perm$perm.p.adj.signif[perm$perm.p.adj < 0.05] = '*'
stat.test = perm %>%
dplyr::select(transcript, group1, group2, perm.p, perm.p.adj, perm.p.adj.signif) %>%
dplyr::rename(
p = perm.p,
p.adj = perm.p.adj,
p.adj.signif = perm.p.adj.signif
)
# Compute y-position using both mean ± se and scaled values
y_pos_df = mydata.long.SE %>%
dplyr::group_by(transcript) %>%
dplyr::summarise(
max_mean_se = max(mean + se, na.rm = TRUE),
min_mean_se = min(mean - se, na.rm = TRUE),
.groups = "drop"
) %>%
dplyr::left_join(
mydata.long %>%
dplyr::group_by(transcript) %>%
dplyr::summarise(
max_scaled = max(scaled, na.rm = TRUE),
min_scaled = min(scaled, na.rm = TRUE),
.groups = "drop"
),
by = "transcript"
) %>%
dplyr::mutate(
plot_max = pmax(max_mean_se, max_scaled, na.rm = TRUE),
plot_min = pmin(min_mean_se, min_scaled, na.rm = TRUE),
plot_range = pmax(plot_max - plot_min, 1e-6)
) %>%
dplyr::select(transcript, plot_max, plot_range)
# Map group names to numeric x positions
pos_map = mydata.long %>%
dplyr::distinct(transcript, Treatment) %>%
dplyr::mutate(xpos = as.numeric(factor(Treatment, levels = levels(mydata.long$Treatment))))
# Attach y-position and x-position info to stat.test
stat.test = stat.test %>%
dplyr::left_join(y_pos_df, by = "transcript") %>%
dplyr::left_join(pos_map %>% dplyr::rename(group1 = Treatment, xmin = xpos), by = c("transcript", "group1")) %>%
dplyr::left_join(pos_map %>% dplyr::rename(group2 = Treatment, xmax = xpos), by = c("transcript", "group2")) %>%
dplyr::mutate(
xmin = dplyr::if_else(is.na(xmin), as.numeric(factor(group1, levels = levels(mydata.long$Treatment))), xmin),
xmax = dplyr::if_else(is.na(xmax), as.numeric(factor(group2, levels = levels(mydata.long$Treatment))), xmax)
)
# Compute distinct y.position per comparison for significant results
stat.test.sig = stat.test %>%
dplyr::filter(!is.na(p.adj) & p.adj <= 0.05) %>%
dplyr::group_by(transcript) %>%
dplyr::arrange(xmin, xmax, .by_group = TRUE) %>%
dplyr::mutate(
base_y = plot_max + 0.03 * plot_range,
inc = 0.05 * plot_range,
y.position = base_y + (dplyr::row_number() - 1) * inc
) %>%
dplyr::ungroup() %>%
dplyr::select(transcript, group1, group2, p.adj, p.adj.signif, y.position, xmin, xmax)
group1 = sapply(y_scales1, function(x) {
s = deparse(x)
s_full = paste(s, collapse = " ")
sub('.*transcript == *"([^"]+)".*', '\\1', s_full)
})
group2 = sapply(y_scales2, function(x) {
s = deparse(x)
s_full = paste(s, collapse = " ")
sub('.*transcript == *"([^"]+)".*', '\\1', s_full)
})
p1 = plot_gene_expression_func(
data.SE = mydata.long.SE,
data.long = mydata.long,
stat.test.sig = stat.test.sig,
transcripts_excl = group2,
facet_cols = 2,
color_values = pal,
plot_title = what,
y_scales = y_scales1,
dodge_width = 0.8
)
print(p1)
p2 = plot_gene_expression_func(
data.SE = mydata.long.SE,
data.long = mydata.long,
stat.test.sig = stat.test.sig,
transcripts_excl = group1,
facet_cols = 6,
color_values = pal,
plot_title = what,
y_scales = y_scales2,
dodge_width = 0.8
)
print(p2)
return(list(plot1 = p1, plot2 = p2, stat.test = stat.test))
}4 Test
4.1 input
## [1] "AMFcolonisation.xlsx" "qPCRmycorrhiza.xlsx" "README.md"
## [4] "Table S1.xlsx" "Table S10.xlsx" "Table S12.xlsx"
## [7] "Table S2.xlsx" "Table S3.xlsx" "Table S4.xlsx"
## [10] "Table S7.xlsx"fn = 'Table S7.xlsx'
PS216 = openxlsx::read.xlsx(xlsxFile = file.path(fp, fn),
sheet = 'Test',
startRow = 1,
colNames = TRUE,
rowNames = FALSE,
detectDates = FALSE,
skipEmptyRows = TRUE,
skipEmptyCols = TRUE,
rows = NULL,
cols = NULL,
check.names = FALSE,
sep.names = ".",
namedRegion = NULL,
na.strings = "NA",
fillMergedCells = FALSE)
data.table::setDT(PS216)
PS216[, Sample.ID := NULL]
PS216[, Genotype := NULL]
PS216$Tissue = as.factor(trimws(PS216$Tissue))
PS216$Treatment = factor(trimws(PS216$Treatment), levels = c("noninoculated", "inoculated"))
PS218 = PS216[grep('PS-216', PS216$Strain, invert = TRUE), ]
PS216 = PS216[grep('PS-218', PS216$Strain, invert = TRUE), ]
PS216[, Strain := NULL]
PS218[, Strain := NULL]4.2 params
4.3 results
shoots.216 = process_data(PS216, groupvars, measurevar, scale_treatment = "noninoculated")$shoots
roots.216 = process_data(PS216, groupvars, measurevar, scale_treatment = "noninoculated")$roots
shoots.218 = process_data(PS218, groupvars, measurevar, scale_treatment = "noninoculated")$shoots
roots.218 = process_data(PS218, groupvars, measurevar, scale_treatment = "noninoculated")$roots
cat("Test", file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")4.3.1 shoots 216
max_per_transcript = shoots.216 %>%
dplyr::group_by(transcript) %>%
dplyr::summarise(max_scaled = max(scaled, na.rm = TRUE), .groups = "drop")
group2_names = setdiff(max_per_transcript$transcript, group1_names)
y_scales1 = build_y_scales_for(group1_names, max_per_transcript)
y_scales2 = build_y_scales_for(group2_names, max_per_transcript)
result = test_and_plot(data_long_raw = shoots.216
, myorder = myorder
, pal = pal
, what = "shoots 216"
, plot_gene_expression_func = plot_gene_expression
, groupvars = c("Treatment", "transcript")
, y_scales1
, y_scales2
)
## #### ####
## Distribution tests
## #### ###### # A tibble: 8 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 8 RBOHD 0.829 0.749 0.804
## 2 8 PR1B 0.00129 0.00125 0.00813
## 3 8 CPI8 0.829 0.749 0.804
## 4 8 CAB 0.829 0.749 0.804
## 5 8 BGLU2 0.829 0.749 0.804
## 6 8 HSP70 0.829 0.749 0.804
## 7 8 PTI5 0.495 0.427 0.264
## 8 8 13-LOX 0.0905 0.108 0.131
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 8 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 RBOHD 1 6 0.181 0.685
## 2 PR1B 1 6 0.478 0.515
## 3 CPI8 1 6 0.725 0.427
## 4 CAB 1 6 0.00927 0.926
## 5 BGLU2 1 6 0.739 0.423
## 6 HSP70 1 6 0.0138 0.910
## 7 PTI5 1 6 1.69 0.242
## 8 13-LOX 1 6 2.61 0.157
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 8 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 RBOHD 1 6 0.181 0.685
## 2 PR1B 1 6 0.478 0.515
## 3 CPI8 1 6 0.725 0.427
## 4 CAB 1 6 0.00927 0.926
## 5 BGLU2 1 6 0.739 0.423
## 6 HSP70 1 6 0.0138 0.910
## 7 PTI5 1 6 1.69 0.242
## 8 13-LOX 1 6 2.61 0.157
## #### ####
## Fligner
## #### ####
## # A tibble: 8 × 5
## # Groups: transcript [8]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 RBOHD 0.0955 0.757 1 Fligner-Killeen test of homogeneity of…
## 2 PR1B 0.0945 0.759 1 Fligner-Killeen test of homogeneity of…
## 3 CPI8 0.120 0.729 1 Fligner-Killeen test of homogeneity of…
## 4 CAB 0.0955 0.757 1 Fligner-Killeen test of homogeneity of…
## 5 BGLU2 0.123 0.726 1 Fligner-Killeen test of homogeneity of…
## 6 HSP70 0.0937 0.760 1 Fligner-Killeen test of homogeneity of…
## 7 PTI5 1.45 0.228 1 Fligner-Killeen test of homogeneity of…
## 8 13-LOX 1.46 0.227 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 8 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.102 RBOHD 4 4 small
## 2 measurement noninoculated inoculated 0.306 PR1B 4 4 moderate
## 3 measurement noninoculated inoculated 0.510 CPI8 4 4 large
## 4 measurement noninoculated inoculated 0.102 CAB 4 4 small
## 5 measurement noninoculated inoculated 0.612 BGLU2 4 4 large
## 6 measurement noninoculated inoculated 0.204 HSP70 4 4 small
## 7 measurement noninoculated inoculated 0.816 PTI5 4 4 large
## 8 measurement noninoculated inoculated 0.714 13-LOX 4 4 large
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 8 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated -0.269 RBOHD 4 4 small
## 2 measurement noninoculated inoculated -0.696 PR1B 4 4 moderate
## 3 measurement noninoculated inoculated -1.21 CPI8 4 4 large
## 4 measurement noninoculated inoculated -0.425 CAB 4 4 small
## 5 measurement noninoculated inoculated -1.44 BGLU2 4 4 large
## 6 measurement noninoculated inoculated -0.389 HSP70 4 4 small
## 7 measurement noninoculated inoculated -4.45 PTI5 4 4 large
## 8 measurement noninoculated inoculated -1.68 13-LOX 4 4 large
res = result$stat.test[, grep("transcript|group1|group2|^p$|p\\.", colnames(result$stat.test))]
print(res[res$p.adj.signif != 'ns', ])## transcript group1 group2 p p.adj p.adj.signif
## 7 PTI5 noninoculated inoculated 0.003937008 0.03149606 *cat("", file = fr, append = TRUE, sep = "\n")
output_text = "permutational t-test"
cat(output_text, file = fr, append = TRUE, sep = "\n")
header = paste(colnames(res), collapse = "\t")
cat(header, file = fr, append = TRUE, sep = "\n")
output_text = apply(as.data.frame(tibble::as_tibble(res)), 1, function(row) paste(row, collapse = "\t"))
cat(output_text, file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
ggsave(
filename = file.path(my_dir, "Test_shoots.216_1.pdf"),
plot = result$plot1,
device = pdf,
path = NULL,
scale = 1,
width = 3,
height = 8,
units = c("in"),
dpi = 900,
limitsize = TRUE,
bg = NULL
)
ggsave(
filename = file.path(my_dir, "Test_shoots.216_2.pdf"),
plot = result$plot2,
device = pdf,
path = NULL,
scale = 1,
width = 9,
height = 8,
units = c("in"),
dpi = 900,
limitsize = TRUE,
bg = NULL
)4.3.2 roots 216
max_per_transcript = roots.216 %>%
dplyr::group_by(transcript) %>%
dplyr::summarise(max_scaled = max(scaled, na.rm = TRUE), .groups = "drop")
group2_names = setdiff(max_per_transcript$transcript, group1_names)
y_scales1 = build_y_scales_for(group1_names, max_per_transcript)
y_scales2 = build_y_scales_for(group2_names, max_per_transcript)
result = test_and_plot(data_long_raw = roots.216
, myorder = myorder
, pal = pal
, what = "roots 216"
, plot_gene_expression_func = plot_gene_expression
, groupvars = c("Treatment", "transcript")
, y_scales1
, y_scales2
)
## #### ####
## Distribution tests
## #### ####
## # A tibble: 8 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 8 RBOHD 0.500 0.532 0.619
## 2 8 PR1B 0.000177 0.000309 0.000230
## 3 8 CPI8 0.0915 0.121 0.258
## 4 8 CAB 0.0543 0.0534 0.0216
## 5 8 BGLU2 0.000188 0.000339 0.000616
## 6 8 HSP70 0.765 0.809 0.956
## 7 8 PTI5 0.0915 0.126 0.258
## 8 8 13-LOX 0.500 0.532 0.580
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 8 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 RBOHD 1 6 1.19 0.317
## 2 PR1B 1 6 1.00 0.356
## 3 CPI8 1 6 7.39 0.0347
## 4 CAB 1 6 3.90 0.0956
## 5 BGLU2 1 6 0.792 0.408
## 6 HSP70 1 6 0.0123 0.915
## 7 PTI5 1 6 9.03 0.0238
## 8 13-LOX 1 6 0.258 0.630
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 8 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 RBOHD 1 6 1.19 0.317
## 2 PR1B 1 6 1.00 0.356
## 3 CPI8 1 6 7.39 0.0347
## 4 CAB 1 6 3.90 0.0956
## 5 BGLU2 1 6 0.792 0.408
## 6 HSP70 1 6 0.0123 0.915
## 7 PTI5 1 6 9.03 0.0238
## 8 13-LOX 1 6 0.258 0.630
## #### ####
## Fligner
## #### ####
## # A tibble: 8 × 5
## # Groups: transcript [8]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 RBOHD 1.47 0.226 1 Fligner-Killeen test of homogeneity of…
## 2 PR1B 4.99 0.0255 1 Fligner-Killeen test of homogeneity of…
## 3 CPI8 2.85 0.0913 1 Fligner-Killeen test of homogeneity of…
## 4 CAB 5.00 0.0253 1 Fligner-Killeen test of homogeneity of…
## 5 BGLU2 0.619 0.431 1 Fligner-Killeen test of homogeneity of…
## 6 HSP70 0.624 0.430 1 Fligner-Killeen test of homogeneity of…
## 7 PTI5 4.96 0.0259 1 Fligner-Killeen test of homogeneity of…
## 8 13-LOX 0.000202 0.989 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 8 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.306 RBOHD 4 4 moderate
## 2 measurement noninoculated inoculated 0.102 PR1B 4 4 small
## 3 measurement noninoculated inoculated 0.408 CPI8 4 4 moderate
## 4 measurement noninoculated inoculated 0.408 CAB 4 4 moderate
## 5 measurement noninoculated inoculated 0.306 BGLU2 4 4 moderate
## 6 measurement noninoculated inoculated 0.816 HSP70 4 4 large
## 7 measurement noninoculated inoculated 0.816 PTI5 4 4 large
## 8 measurement noninoculated inoculated 0.816 13-LOX 4 4 large
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 8 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.318 RBOHD 4 4 small
## 2 measurement noninoculated inoculated -0.653 PR1B 4 4 moderate
## 3 measurement noninoculated inoculated -1.26 CPI8 4 4 large
## 4 measurement noninoculated inoculated -1.60 CAB 4 4 large
## 5 measurement noninoculated inoculated -0.499 BGLU2 4 4 small
## 6 measurement noninoculated inoculated -2.82 HSP70 4 4 large
## 7 measurement noninoculated inoculated -2.43 PTI5 4 4 large
## 8 measurement noninoculated inoculated -3.09 13-LOX 4 4 large
res = result$stat.test[, grep("transcript|group1|group2|^p$|p\\.", colnames(result$stat.test))]
print(res[res$p.adj.signif != 'ns', ])## transcript group1 group2 p p.adj p.adj.signif
## 6 HSP70 noninoculated inoculated 0.003937008 0.01049869 *
## 7 PTI5 noninoculated inoculated 0.003937008 0.01049869 *
## 8 13-LOX noninoculated inoculated 0.003937008 0.01049869 *cat("", file = fr, append = TRUE, sep = "\n")
output_text = "permutational t-test"
cat(output_text, file = fr, append = TRUE, sep = "\n")
header = paste(colnames(res), collapse = "\t")
cat(header, file = fr, append = TRUE, sep = "\n")
output_text = apply(as.data.frame(tibble::as_tibble(res)), 1, function(row) paste(row, collapse = "\t"))
cat(output_text, file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
ggsave(
filename = file.path(my_dir, "Test_roots.216_1.pdf"),
plot = result$plot1,
device = pdf,
path = NULL,
scale = 1,
width = 3,
height = 8,
units = c("in"),
dpi = 900,
limitsize = TRUE,
bg = NULL
)
ggsave(
filename = file.path(my_dir, "Test_roots.216_2.pdf"),
plot = result$plot2,
device = pdf,
path = NULL,
scale = 1,
width = 9,
height = 8,
units = c("in"),
dpi = 900,
limitsize = TRUE,
bg = NULL
)4.3.3 shoots 218
max_per_transcript = shoots.218 %>%
dplyr::group_by(transcript) %>%
dplyr::summarise(max_scaled = max(scaled, na.rm = TRUE), .groups = "drop")
group2_names = setdiff(max_per_transcript$transcript, group1_names)
y_scales1 = build_y_scales_for(group1_names, max_per_transcript)
y_scales2 = build_y_scales_for(group2_names, max_per_transcript)
result = test_and_plot(data_long_raw = shoots.218
, myorder = myorder
, pal = pal
, what = "shoots 218"
, plot_gene_expression_func = plot_gene_expression
, groupvars = c("Treatment", "transcript")
, y_scales1
, y_scales2
)
## #### ####
## Distribution tests
## #### ####
## # A tibble: 8 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 8 RBOHD 0.658 0.495 0.634
## 2 8 PR1B 0.000437 0.000450 0.00552
## 3 8 CPI8 0.0280 0.0266 0.0614
## 4 8 CAB 0.658 0.489 0.604
## 5 8 BGLU2 0.522 0.489 0.604
## 6 8 HSP70 0.658 0.505 0.604
## 7 8 PTI5 0.522 0.489 0.604
## 8 8 13-LOX 0.522 0.489 0.604
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 8 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 RBOHD 1 6 0.255 0.631
## 2 PR1B 1 6 0.585 0.473
## 3 CPI8 1 6 3.74 0.101
## 4 CAB 1 6 0.0364 0.855
## 5 BGLU2 1 6 0.00507 0.946
## 6 HSP70 1 6 0.000737 0.979
## 7 PTI5 1 6 0.0321 0.864
## 8 13-LOX 1 6 1.91 0.216
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 8 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 RBOHD 1 6 0.255 0.631
## 2 PR1B 1 6 0.585 0.473
## 3 CPI8 1 6 3.74 0.101
## 4 CAB 1 6 0.0364 0.855
## 5 BGLU2 1 6 0.00507 0.946
## 6 HSP70 1 6 0.000737 0.979
## 7 PTI5 1 6 0.0321 0.864
## 8 13-LOX 1 6 1.91 0.216
## #### ####
## Fligner
## #### ####
## # A tibble: 8 × 5
## # Groups: transcript [8]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 RBOHD 0.715 0.398 1 Fligner-Killeen test of homogeneity of…
## 2 PR1B 0.303 0.582 1 Fligner-Killeen test of homogeneity of…
## 3 CPI8 5.03 0.0249 1 Fligner-Killeen test of homogeneity of…
## 4 CAB 0.0947 0.758 1 Fligner-Killeen test of homogeneity of…
## 5 BGLU2 0.623 0.430 1 Fligner-Killeen test of homogeneity of…
## 6 HSP70 0.0955 0.757 1 Fligner-Killeen test of homogeneity of…
## 7 PTI5 0.0955 0.757 1 Fligner-Killeen test of homogeneity of…
## 8 13-LOX 0.709 0.400 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 8 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.408 RBOHD 4 4 moderate
## 2 measurement noninoculated inoculated 0.204 PR1B 4 4 small
## 3 measurement noninoculated inoculated 0 CPI8 4 4 small
## 4 measurement noninoculated inoculated 0.306 CAB 4 4 moderate
## 5 measurement noninoculated inoculated 0.612 BGLU2 4 4 large
## 6 measurement noninoculated inoculated 0.612 HSP70 4 4 large
## 7 measurement noninoculated inoculated 0.816 PTI5 4 4 large
## 8 measurement noninoculated inoculated 0.612 13-LOX 4 4 large
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 8 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated -1.23 RBOHD 4 4 large
## 2 measurement noninoculated inoculated -0.702 PR1B 4 4 moderate
## 3 measurement noninoculated inoculated 0.436 CPI8 4 4 small
## 4 measurement noninoculated inoculated -0.615 CAB 4 4 moderate
## 5 measurement noninoculated inoculated -1.79 BGLU2 4 4 large
## 6 measurement noninoculated inoculated -1.77 HSP70 4 4 large
## 7 measurement noninoculated inoculated -4.62 PTI5 4 4 large
## 8 measurement noninoculated inoculated -1.37 13-LOX 4 4 large
res = result$stat.test[, grep("transcript|group1|group2|^p$|p\\.", colnames(result$stat.test))]
print(res[res$p.adj.signif != 'ns', ])## [1] transcript group1 group2 p p.adj
## [6] p.adj.signif
## <0 rows> (or 0-length row.names)cat("", file = fr, append = TRUE, sep = "\n")
output_text = "permutational t-test"
cat(output_text, file = fr, append = TRUE, sep = "\n")
header = paste(colnames(res), collapse = "\t")
cat(header, file = fr, append = TRUE, sep = "\n")
output_text = apply(as.data.frame(tibble::as_tibble(res)), 1, function(row) paste(row, collapse = "\t"))
cat(output_text, file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
ggsave(
filename = file.path(my_dir, "Test_shoots.218_1.pdf"),
plot = result$plot1,
device = pdf,
path = NULL,
scale = 1,
width = 3,
height = 8,
units = c("in"),
dpi = 900,
limitsize = TRUE,
bg = NULL
)
ggsave(
filename = file.path(my_dir, "Test_shoots.218_2.pdf"),
plot = result$plot2,
device = pdf,
path = NULL,
scale = 1,
width = 9,
height = 8,
units = c("in"),
dpi = 900,
limitsize = TRUE,
bg = NULL
)4.3.4 roots 218
max_per_transcript = roots.218 %>%
dplyr::group_by(transcript) %>%
dplyr::summarise(max_scaled = max(scaled, na.rm = TRUE), .groups = "drop")
group2_names = setdiff(max_per_transcript$transcript, group1_names)
y_scales1 = build_y_scales_for(group1_names, max_per_transcript)
y_scales2 = build_y_scales_for(group2_names, max_per_transcript)
result = test_and_plot(data_long_raw = roots.218
, myorder = myorder
, pal = pal
, what = "roots 218"
, plot_gene_expression_func = plot_gene_expression
, groupvars = c("Treatment", "transcript")
, y_scales1
, y_scales2
)
## #### ####
## Distribution tests
## #### ####
## # A tibble: 8 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 8 RBOHD 0.271 0.322 0.573
## 2 8 PR1B 0.00339 0.00629 0.00158
## 3 8 CPI8 0.0799 0.115 0.296
## 4 8 CAB 0.0650 0.0596 0.0195
## 5 8 BGLU2 0.320 0.348 0.573
## 6 8 HSP70 0.320 0.335 0.539
## 7 8 PTI5 0.107 0.121 0.296
## 8 8 13-LOX 0.327 0.378 0.573
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 8 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 RBOHD 1 6 8.01 0.0300
## 2 PR1B 1 6 0.668 0.445
## 3 CPI8 1 6 2.55 0.161
## 4 CAB 1 6 6.52 0.0433
## 5 BGLU2 1 6 1.96 0.211
## 6 HSP70 1 6 0.0221 0.887
## 7 PTI5 1 6 24.4 0.00260
## 8 13-LOX 1 6 4.36 0.0818
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 8 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 RBOHD 1 6 8.01 0.0300
## 2 PR1B 1 6 0.668 0.445
## 3 CPI8 1 6 2.55 0.161
## 4 CAB 1 6 6.52 0.0433
## 5 BGLU2 1 6 1.96 0.211
## 6 HSP70 1 6 0.0221 0.887
## 7 PTI5 1 6 24.4 0.00260
## 8 13-LOX 1 6 4.36 0.0818
## #### ####
## Fligner
## #### ####
## # A tibble: 8 × 5
## # Groups: transcript [8]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 RBOHD 5.03 0.0249 1 Fligner-Killeen test of homogeneity of…
## 2 PR1B 0.123 0.726 1 Fligner-Killeen test of homogeneity of…
## 3 CPI8 2.85 0.0913 1 Fligner-Killeen test of homogeneity of…
## 4 CAB 4.96 0.0259 1 Fligner-Killeen test of homogeneity of…
## 5 BGLU2 1.75 0.186 1 Fligner-Killeen test of homogeneity of…
## 6 HSP70 0.000272 0.987 1 Fligner-Killeen test of homogeneity of…
## 7 PTI5 5.00 0.0253 1 Fligner-Killeen test of homogeneity of…
## 8 13-LOX 1.76 0.185 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 8 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.306 RBOHD 4 4 moderate
## 2 measurement noninoculated inoculated 0.408 PR1B 4 4 moderate
## 3 measurement noninoculated inoculated 0.510 CPI8 4 4 large
## 4 measurement noninoculated inoculated 0.816 CAB 4 4 large
## 5 measurement noninoculated inoculated 0.408 BGLU2 4 4 moderate
## 6 measurement noninoculated inoculated 0.204 HSP70 4 4 small
## 7 measurement noninoculated inoculated 0.816 PTI5 4 4 large
## 8 measurement noninoculated inoculated 0.816 13-LOX 4 4 large
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 8 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated -0.928 RBOHD 4 4 large
## 2 measurement noninoculated inoculated -0.769 PR1B 4 4 moderate
## 3 measurement noninoculated inoculated -1.44 CPI8 4 4 large
## 4 measurement noninoculated inoculated -1.94 CAB 4 4 large
## 5 measurement noninoculated inoculated 0.906 BGLU2 4 4 large
## 6 measurement noninoculated inoculated -0.560 HSP70 4 4 moderate
## 7 measurement noninoculated inoculated -2.12 PTI5 4 4 large
## 8 measurement noninoculated inoculated -2.81 13-LOX 4 4 large
res = result$stat.test[, grep("transcript|group1|group2|^p$|p\\.", colnames(result$stat.test))]
print(res[res$p.adj.signif != 'ns', ])## transcript group1 group2 p p.adj p.adj.signif
## 4 CAB noninoculated inoculated 0.003937008 0.01049869 *
## 7 PTI5 noninoculated inoculated 0.003937008 0.01049869 *
## 8 13-LOX noninoculated inoculated 0.003937008 0.01049869 *cat("", file = fr, append = TRUE, sep = "\n")
output_text = "permutational t-test"
cat(output_text, file = fr, append = TRUE, sep = "\n")
header = paste(colnames(res), collapse = "\t")
cat(header, file = fr, append = TRUE, sep = "\n")
output_text = apply(as.data.frame(tibble::as_tibble(res)), 1, function(row) paste(row, collapse = "\t"))
cat(output_text, file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
ggsave(
filename = file.path(my_dir, "Test_roots.218_1.pdf"),
plot = result$plot1,
device = pdf,
path = NULL,
scale = 1,
width = 3,
height = 8,
units = c("in"),
dpi = 900,
limitsize = TRUE,
bg = NULL
)
ggsave(
filename = file.path(my_dir, "Test_roots.218_2.pdf"),
plot = result$plot2,
device = pdf,
path = NULL,
scale = 1,
width = 9,
height = 8,
units = c("in"),
dpi = 900,
limitsize = TRUE,
bg = NULL
)5 Exp1
5.1 input
## [1] "AMFcolonisation.xlsx" "qPCRmycorrhiza.xlsx" "README.md"
## [4] "Table S1.xlsx" "Table S10.xlsx" "Table S12.xlsx"
## [7] "Table S2.xlsx" "Table S3.xlsx" "Table S4.xlsx"
## [10] "Table S7.xlsx"fn = 'Table S7.xlsx'
PS218 = openxlsx::read.xlsx(xlsxFile = file.path(fp, fn),
sheet = 'Exp1',
startRow = 1,
colNames = TRUE,
rowNames = FALSE,
detectDates = FALSE,
skipEmptyRows = TRUE,
skipEmptyCols = TRUE,
rows = NULL,
cols = NULL,
check.names = FALSE,
sep.names = ".",
namedRegion = NULL,
na.strings = "NA",
fillMergedCells = FALSE)
data.table::setDT(PS218)
PS218[, Sample.ID := NULL]
PS218[, Genotype := NULL]
PS218[, Strain := NULL]
PS218$Tissue = as.factor(trimws(PS218$Tissue))
PS218$Treatment = factor(trimws(PS218$Treatment), levels = c("noninoculated", "inoculated"))
table(PS218$Time)##
## 2 h 4 h 6 h
## 16 16 165.2 params
5.3 results
cat("EXP1", file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
run_analysis_for_time <- function(data
, time_point
, myorder
, pal
, groupvars
, measurevar
, my_dir
, plot_gene_expression
, test_and_plot) {
temp = data[data$Time == time_point, ]
temp[, Time := NULL]
shoots = process_data(temp, groupvars, measurevar, scale_treatment = "noninoculated")$shoots
roots = process_data(temp, groupvars, measurevar, scale_treatment = "noninoculated")$roots
max_per_transcript = shoots %>%
dplyr::group_by(transcript) %>%
dplyr::summarise(max_scaled = max(scaled, na.rm = TRUE), .groups = "drop")
group2_names = setdiff(max_per_transcript$transcript, group1_names)
y_scales1 = build_y_scales_for(group1_names, max_per_transcript)
y_scales2 = build_y_scales_for(group2_names, max_per_transcript)
max_per_transcript = roots %>%
dplyr::group_by(transcript) %>%
dplyr::summarise(max_scaled = max(scaled, na.rm = TRUE), .groups = "drop")
y_scales3 = build_y_scales_for(group1_names, max_per_transcript)
y_scales4 = build_y_scales_for(group2_names, max_per_transcript)
result_shoots = test_and_plot(data_long_raw = shoots,
myorder = myorder,
pal = pal,
what = paste("shoots", time_point),
plot_gene_expression_func = plot_gene_expression,
groupvars = groupvars,
y_scales1 = y_scales1,
y_scales2 = y_scales2)
res = result_shoots$stat.test[, grep("transcript|group1|group2|^p$|p\\.", colnames(result_shoots$stat.test))]
print(res[res$p.adj.signif != 'ns', ])
cat("", file = fr, append = TRUE, sep = "\n")
output_text = "permutational t-test"
cat(output_text, file = fr, append = TRUE, sep = "\n")
header = paste(colnames(res), collapse = "\t")
cat(header, file = fr, append = TRUE, sep = "\n")
output_text = apply(as.data.frame(tibble::as_tibble(res)), 1, function(row) paste(row, collapse = "\t"))
cat(output_text, file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
ggsave(filename = file.path(my_dir, paste0("Exp1_shoots.", gsub(" ", "", time_point), "_1.pdf")),
plot = result_shoots$plot1, device = pdf, width = 3, height = 8, units = "in", dpi = 900)
ggsave(filename = file.path(my_dir, paste0("Exp1_shoots.", gsub(" ", "", time_point), "_2.pdf")),
plot = result_shoots$plot2, device = pdf, width = 3, height = 8, units = "in", dpi = 900)
result_roots = test_and_plot(data_long_raw = roots,
myorder = myorder,
pal = pal,
what = paste("roots", time_point),
plot_gene_expression_func = plot_gene_expression,
groupvars = groupvars,
y_scales1 = y_scales3,
y_scales2 = y_scales4)
res = result_roots$stat.test[, grep("transcript|group1|group2|^p$|p\\.", colnames(result_roots$stat.test))]
print(res[res$p.adj.signif != 'ns', ])
cat("", file = fr, append = TRUE, sep = "\n")
output_text = "permutational t-test"
cat(output_text, file = fr, append = TRUE, sep = "\n")
header = paste(colnames(res), collapse = "\t")
cat(header, file = fr, append = TRUE, sep = "\n")
output_text = apply(as.data.frame(tibble::as_tibble(res)), 1, function(row) paste(row, collapse = "\t"))
cat(output_text, file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
ggsave(filename = file.path(my_dir, paste0("Exp1_roots.", gsub(" ", ".", time_point), "_1.pdf")),
plot = result_roots$plot1, device = pdf, width = 3, height = 8, units = "in", dpi = 900)
ggsave(filename = file.path(my_dir, paste0("Exp1_roots.", gsub(" ", ".", time_point), "_2.pdf")),
plot = result_roots$plot2, device = pdf, width = 3, height = 8, units = "in", dpi = 900)
list(shoots = result_shoots, roots = result_roots)
}
results_2h = run_analysis_for_time(data = PS218
, time_point = "2 h"
, myorder
, pal
, groupvars
, measurevar
, my_dir
, plot_gene_expression
, test_and_plot)
## #### ####
## Distribution tests
## #### ####
## # A tibble: 4 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 8 BGLU2 0.773 0.687 0.574
## 2 8 HSP70 0.773 0.687 0.574
## 3 8 PTI5 0.773 0.687 0.631
## 4 8 13-LOX 0.773 0.696 0.689
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 1.91 0.216
## 2 HSP70 1 6 0.0663 0.805
## 3 PTI5 1 6 0.156 0.707
## 4 13-LOX 1 6 0.0281 0.872
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 1.91 0.216
## 2 HSP70 1 6 0.0663 0.805
## 3 PTI5 1 6 0.156 0.707
## 4 13-LOX 1 6 0.0281 0.872
## #### ####
## Fligner
## #### ####
## # A tibble: 4 × 5
## # Groups: transcript [4]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 BGLU2 1.76 0.185 1 Fligner-Killeen test of homogeneity of…
## 2 HSP70 0.619 0.431 1 Fligner-Killeen test of homogeneity of…
## 3 PTI5 0.303 0.582 1 Fligner-Killeen test of homogeneity of…
## 4 13-LOX 0.120 0.729 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.510 BGLU2 4 4 large
## 2 measurement noninoculated inoculated 0.102 HSP70 4 4 small
## 3 measurement noninoculated inoculated 0.612 PTI5 4 4 large
## 4 measurement noninoculated inoculated 0 13-LOX 4 4 small
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 1.13 BGLU2 4 4 large
## 2 measurement noninoculated inoculated -0.176 HSP70 4 4 negligible
## 3 measurement noninoculated inoculated -1.19 PTI5 4 4 large
## 4 measurement noninoculated inoculated 0.0435 13-LOX 4 4 negligible
## [1] transcript group1 group2 p p.adj
## [6] p.adj.signif
## <0 rows> (or 0-length row.names)
## #### ####
## Distribution tests
## #### ####
## # A tibble: 4 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 7 BGLU2 0.252 NA 0.717
## 2 8 HSP70 0.0417 0.0347 0.148
## 3 8 PTI5 0.0538 0.0606 0.356
## 4 8 13-LOX 0.0417 0.0347 0.144
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 5 2.63 0.166
## 2 HSP70 1 6 38.9 0.000785
## 3 PTI5 1 6 4.06 0.0906
## 4 13-LOX 1 6 0.370 0.566
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 5 2.63 0.166
## 2 HSP70 1 6 38.9 0.000785
## 3 PTI5 1 6 4.06 0.0906
## 4 13-LOX 1 6 0.370 0.566
## #### ####
## Fligner
## #### ####
## # A tibble: 4 × 5
## # Groups: transcript [4]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 BGLU2 0.923 0.337 1 Fligner-Killeen test of homogeneity of…
## 2 HSP70 4.99 0.0255 1 Fligner-Killeen test of homogeneity of…
## 3 PTI5 5.00 0.0253 1 Fligner-Killeen test of homogeneity of…
## 4 13-LOX 0.0947 0.758 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.267 BGLU2 4 3 small
## 2 measurement noninoculated inoculated 0.204 HSP70 4 4 small
## 3 measurement noninoculated inoculated 0.816 PTI5 4 4 large
## 4 measurement noninoculated inoculated 0.816 13-LOX 4 4 large
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated -0.789 BGLU2 4 3 moderate
## 2 measurement noninoculated inoculated -0.823 HSP70 4 4 large
## 3 measurement noninoculated inoculated -1.99 PTI5 4 4 large
## 4 measurement noninoculated inoculated -11.7 13-LOX 4 4 large
## transcript group1 group2 p p.adj p.adj.signif
## 3 PTI5 noninoculated inoculated 0.003937008 0.007874016 *
## 4 13-LOX noninoculated inoculated 0.003937008 0.007874016 *
results_4h = run_analysis_for_time(data = PS218
, time_point = "4 h"
, myorder
, pal
, groupvars
, measurevar
, my_dir
, plot_gene_expression
, test_and_plot)
## #### ####
## Distribution tests
## #### ####
## # A tibble: 4 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 8 BGLU2 0.424 0.492 0.703
## 2 8 HSP70 0.424 0.492 0.723
## 3 7 PTI5 0.0396 NA 0.0196
## 4 8 13-LOX 0.424 0.492 0.723
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 7.26 0.0358
## 2 HSP70 1 6 2.06 0.201
## 3 PTI5 1 5 3.62 0.115
## 4 13-LOX 1 6 0.0744 0.794
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 7.26 0.0358
## 2 HSP70 1 6 2.06 0.201
## 3 PTI5 1 5 3.62 0.115
## 4 13-LOX 1 6 0.0744 0.794
## #### ####
## Fligner
## #### ####
## # A tibble: 4 × 5
## # Groups: transcript [4]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 BGLU2 2.84 0.0920 1 Fligner-Killeen test of homogeneity of…
## 2 HSP70 2.87 0.0903 1 Fligner-Killeen test of homogeneity of…
## 3 PTI5 0.935 0.334 1 Fligner-Killeen test of homogeneity of…
## 4 13-LOX 0.0945 0.759 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.408 BGLU2 4 4 moderate
## 2 measurement noninoculated inoculated 0.204 HSP70 4 4 small
## 3 measurement noninoculated inoculated 0.668 PTI5 4 3 large
## 4 measurement noninoculated inoculated 0.510 13-LOX 4 4 large
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 1.14 BGLU2 4 4 large
## 2 measurement noninoculated inoculated -0.0196 HSP70 4 4 negligible
## 3 measurement noninoculated inoculated -1.41 PTI5 4 3 large
## 4 measurement noninoculated inoculated -0.802 13-LOX 4 4 large
## [1] transcript group1 group2 p p.adj
## [6] p.adj.signif
## <0 rows> (or 0-length row.names)
## #### ####
## Distribution tests
## #### ####
## # A tibble: 4 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 8 BGLU2 0.343 0.350 0.455
## 2 8 HSP70 0.343 0.350 0.543
## 3 7 PTI5 0.105 NA 0.252
## 4 8 13-LOX 0.343 0.350 0.453
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 0.0398 0.849
## 2 HSP70 1 6 0.494 0.508
## 3 PTI5 1 5 13.7 0.0140
## 4 13-LOX 1 6 2.52 0.164
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 0.0398 0.849
## 2 HSP70 1 6 0.494 0.508
## 3 PTI5 1 5 13.7 0.0140
## 4 13-LOX 1 6 2.52 0.164
## #### ####
## Fligner
## #### ####
## # A tibble: 4 × 5
## # Groups: transcript [4]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 BGLU2 0.0947 0.758 1 Fligner-Killeen test of homogeneity of…
## 2 HSP70 0.715 0.398 1 Fligner-Killeen test of homogeneity of…
## 3 PTI5 4.01 0.0453 1 Fligner-Killeen test of homogeneity of…
## 4 13-LOX 2.85 0.0913 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.204 BGLU2 4 4 small
## 2 measurement noninoculated inoculated 0.714 HSP70 4 4 large
## 3 measurement noninoculated inoculated 0.802 PTI5 3 4 large
## 4 measurement noninoculated inoculated 0.816 13-LOX 4 4 large
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.101 BGLU2 4 4 negligible
## 2 measurement noninoculated inoculated -2.66 HSP70 4 4 large
## 3 measurement noninoculated inoculated -1.62 PTI5 3 4 large
## 4 measurement noninoculated inoculated -2.65 13-LOX 4 4 large
## transcript group1 group2 p p.adj p.adj.signif
## 2 HSP70 noninoculated inoculated 0.035433071 0.04724409 *
## 3 PTI5 noninoculated inoculated 0.023622047 0.04724409 *
## 4 13-LOX noninoculated inoculated 0.003937008 0.01574803 *
results_6h = run_analysis_for_time(data = PS218
, time_point = "6 h"
, myorder
, pal
, groupvars
, measurevar
, my_dir
, plot_gene_expression
, test_and_plot)
## #### ####
## Distribution tests
## #### ####
## # A tibble: 4 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 8 BGLU2 0.410 0.420 0.585
## 2 8 HSP70 0.932 0.825 0.780
## 3 8 PTI5 0.479 0.696 0.780
## 4 8 13-LOX 0.932 0.825 0.780
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 9.33 0.0224
## 2 HSP70 1 6 1.76 0.233
## 3 PTI5 1 6 0.632 0.457
## 4 13-LOX 1 6 0.157 0.705
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 9.33 0.0224
## 2 HSP70 1 6 1.76 0.233
## 3 PTI5 1 6 0.632 0.457
## 4 13-LOX 1 6 0.157 0.705
## #### ####
## Fligner
## #### ####
## # A tibble: 4 × 5
## # Groups: transcript [4]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 BGLU2 4.99 0.0255 1 Fligner-Killeen test of homogeneity of…
## 2 HSP70 1.47 0.226 1 Fligner-Killeen test of homogeneity of…
## 3 PTI5 1.76 0.185 1 Fligner-Killeen test of homogeneity of…
## 4 13-LOX 0.0937 0.760 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.204 BGLU2 4 4 small
## 2 measurement noninoculated inoculated 0.204 HSP70 4 4 small
## 3 measurement noninoculated inoculated 0.714 PTI5 4 4 large
## 4 measurement noninoculated inoculated 0.510 13-LOX 4 4 large
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated -0.783 BGLU2 4 4 moderate
## 2 measurement noninoculated inoculated -0.167 HSP70 4 4 negligible
## 3 measurement noninoculated inoculated -2.26 PTI5 4 4 large
## 4 measurement noninoculated inoculated -1.01 13-LOX 4 4 large
## [1] transcript group1 group2 p p.adj
## [6] p.adj.signif
## <0 rows> (or 0-length row.names)
## #### ####
## Distribution tests
## #### ####
## # A tibble: 4 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 8 BGLU2 0.163 0.184 0.127
## 2 8 HSP70 0.514 0.450 0.330
## 3 8 PTI5 0.163 0.184 0.127
## 4 8 13-LOX 0.163 0.185 0.330
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 1.48 0.269
## 2 HSP70 1 6 0.0375 0.853
## 3 PTI5 1 6 1.80 0.228
## 4 13-LOX 1 6 6.91 0.0391
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 1.48 0.269
## 2 HSP70 1 6 0.0375 0.853
## 3 PTI5 1 6 1.80 0.228
## 4 13-LOX 1 6 6.91 0.0391
## #### ####
## Fligner
## #### ####
## # A tibble: 4 × 5
## # Groups: transcript [4]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 BGLU2 0.123 0.726 1 Fligner-Killeen test of homogeneity of…
## 2 HSP70 0.0947 0.758 1 Fligner-Killeen test of homogeneity of…
## 3 PTI5 1.47 0.226 1 Fligner-Killeen test of homogeneity of…
## 4 13-LOX 5.00 0.0253 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.714 BGLU2 4 4 large
## 2 measurement noninoculated inoculated 0.306 HSP70 4 4 moderate
## 3 measurement noninoculated inoculated 0.816 PTI5 4 4 large
## 4 measurement noninoculated inoculated 0.816 13-LOX 4 4 large
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated -2.10 BGLU2 4 4 large
## 2 measurement noninoculated inoculated -0.886 HSP70 4 4 large
## 3 measurement noninoculated inoculated -4.39 PTI5 4 4 large
## 4 measurement noninoculated inoculated -2.64 13-LOX 4 4 large
## transcript group1 group2 p p.adj p.adj.signif
## 1 BGLU2 noninoculated inoculated 0.023622047 0.031496063 *
## 3 PTI5 noninoculated inoculated 0.003937008 0.007874016 *
## 4 13-LOX noninoculated inoculated 0.003937008 0.007874016 *
6 Exp2
6.1 input
## [1] "AMFcolonisation.xlsx" "qPCRmycorrhiza.xlsx" "README.md"
## [4] "Table S1.xlsx" "Table S10.xlsx" "Table S12.xlsx"
## [7] "Table S2.xlsx" "Table S3.xlsx" "Table S4.xlsx"
## [10] "Table S7.xlsx"fn = 'Table S7.xlsx'
PS218 = openxlsx::read.xlsx(xlsxFile = file.path(fp, fn),
sheet = 'Exp2',
startRow = 1,
colNames = TRUE,
rowNames = FALSE,
detectDates = FALSE,
skipEmptyRows = TRUE,
skipEmptyCols = TRUE,
rows = NULL,
cols = NULL,
check.names = FALSE,
sep.names = ".",
namedRegion = NULL,
na.strings = "NA",
fillMergedCells = FALSE)
data.table::setDT(PS218)
PS218[, Sample.ID := NULL]
PS218[, Genotype := NULL]
PS218[, Strain := NULL]
PS218$Tissue = as.factor(trimws(PS218$Tissue))
PS218$Treatment = factor(trimws(PS218$Treatment), levels = c("noninoculated", "inoculated"))
table(PS218$Time)##
## 1 h 1 min 15 min 2 h 30 min
## 8 8 8 8 86.2 params
6.3 results
run_analysis_for_time <- function(data
, time_point
, myorder
, pal
, groupvars
, measurevar
, y_scales1
, y_scales2
, my_dir
, plot_gene_expression
, test_and_plot
) {
temp = data[data$Time == time_point, ]
temp[, Time := NULL]
roots = process_data(temp, groupvars, measurevar, scale_treatment = "noninoculated")$roots
max_per_transcript = roots %>%
dplyr::group_by(transcript) %>%
dplyr::summarise(max_scaled = max(scaled, na.rm = TRUE), .groups = "drop")
group2_names = setdiff(max_per_transcript$transcript, group1_names)
y_scales1 = build_y_scales_for(group1_names, max_per_transcript)
y_scales2 = build_y_scales_for(group2_names, max_per_transcript)
result_roots = test_and_plot(data_long_raw = roots,
myorder = myorder,
pal = pal,
what = paste("roots", time_point),
plot_gene_expression_func = plot_gene_expression,
groupvars = groupvars,
y_scales1,
y_scales2)
res = result_roots$stat.test[, grep("transcript|group1|group2|^p$|p\\.", colnames(result_roots$stat.test))]
print(res[res$p.adj.signif != 'ns', ])
cat("", file = fr, append = TRUE, sep = "\n")
output_text = "permutational t-test"
cat(output_text, file = fr, append = TRUE, sep = "\n")
header = paste(colnames(res), collapse = "\t")
cat(header, file = fr, append = TRUE, sep = "\n")
output_text = apply(as.data.frame(tibble::as_tibble(res)), 1, function(row) paste(row, collapse = "\t"))
cat(output_text, file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
cat("", file = fr, append = TRUE, sep = "\n")
ggsave(filename = file.path(my_dir, paste0("Exp2_roots.", gsub(" ", ".", time_point), "_1.pdf")),
plot = result_roots$plot1, device = pdf, width = 3, height = 8, units = "in", dpi = 900)
ggsave(filename = file.path(my_dir, paste0("Exp2_roots.", gsub(" ", ".", time_point), "_2.pdf")),
plot = result_roots$plot2, device = pdf, width = 3, height = 8, units = "in", dpi = 900)
list(roots = result_roots)
}
results_1min = run_analysis_for_time(data = PS218
, time_point = "1 min"
, myorder
, pal
, groupvars
, measurevar
, y_scales1
, y_scales2
, my_dir
, plot_gene_expression
, test_and_plot
)
## #### ####
## Distribution tests
## #### ####
## # A tibble: 4 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 8 BGLU2 0.556 0.328 0.343
## 2 8 HSP70 0.774 0.511 0.400
## 3 8 PTI5 0.575 0.328 0.400
## 4 8 13-LOX 0.556 0.328 0.400
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 1.63 0.249
## 2 HSP70 1 6 0.530 0.494
## 3 PTI5 1 6 0.000531 0.982
## 4 13-LOX 1 6 0.652 0.450
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 1.63 0.249
## 2 HSP70 1 6 0.530 0.494
## 3 PTI5 1 6 0.000531 0.982
## 4 13-LOX 1 6 0.652 0.450
## #### ####
## Fligner
## #### ####
## # A tibble: 4 × 5
## # Groups: transcript [4]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 BGLU2 0.123 0.726 1 Fligner-Killeen test of homogeneity of…
## 2 HSP70 0.710 0.399 1 Fligner-Killeen test of homogeneity of…
## 3 PTI5 0.0945 0.759 1 Fligner-Killeen test of homogeneity of…
## 4 13-LOX 0.304 0.581 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.612 BGLU2 4 4 large
## 2 measurement noninoculated inoculated 0.204 HSP70 4 4 small
## 3 measurement noninoculated inoculated 0.102 PTI5 4 4 small
## 4 measurement noninoculated inoculated 0.204 13-LOX 4 4 small
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated -1.05 BGLU2 4 4 large
## 2 measurement noninoculated inoculated -0.798 HSP70 4 4 moderate
## 3 measurement noninoculated inoculated 0.430 PTI5 4 4 small
## 4 measurement noninoculated inoculated 0.137 13-LOX 4 4 negligible
## [1] transcript group1 group2 p p.adj
## [6] p.adj.signif
## <0 rows> (or 0-length row.names)
results_15min = run_analysis_for_time(data = PS218
, time_point = "15 min"
, myorder
, pal
, groupvars
, measurevar
, y_scales1
, y_scales2
, my_dir
, plot_gene_expression
, test_and_plot
)
## #### ####
## Distribution tests
## #### ####
## # A tibble: 4 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 8 BGLU2 0.480 0.666 0.948
## 2 8 HSP70 0.920 0.910 0.948
## 3 8 PTI5 0.920 0.910 0.948
## 4 8 13-LOX 0.920 0.910 0.948
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 2.83 0.143
## 2 HSP70 1 6 0.795 0.407
## 3 PTI5 1 6 0.890 0.382
## 4 13-LOX 1 6 3.67 0.104
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 2.83 0.143
## 2 HSP70 1 6 0.795 0.407
## 3 PTI5 1 6 0.890 0.382
## 4 13-LOX 1 6 3.67 0.104
## #### ####
## Fligner
## #### ####
## # A tibble: 4 × 5
## # Groups: transcript [4]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 BGLU2 2.84 0.0920 1 Fligner-Killeen test of homogeneity of…
## 2 HSP70 1.47 0.226 1 Fligner-Killeen test of homogeneity of…
## 3 PTI5 1.46 0.226 1 Fligner-Killeen test of homogeneity of…
## 4 13-LOX 1.46 0.227 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.306 BGLU2 4 4 moderate
## 2 measurement noninoculated inoculated 0.816 HSP70 4 4 large
## 3 measurement noninoculated inoculated 0.408 PTI5 4 4 moderate
## 4 measurement noninoculated inoculated 0.204 13-LOX 4 4 small
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.662 BGLU2 4 4 moderate
## 2 measurement noninoculated inoculated 2.07 HSP70 4 4 large
## 3 measurement noninoculated inoculated -0.883 PTI5 4 4 large
## 4 measurement noninoculated inoculated 0.623 13-LOX 4 4 moderate
## transcript group1 group2 p p.adj p.adj.signif
## 2 HSP70 noninoculated inoculated 0.003937008 0.01574803 *
results_30min = run_analysis_for_time(data = PS218
, time_point = "30 min"
, myorder
, pal
, groupvars
, measurevar
, y_scales1
, y_scales2
, my_dir
, plot_gene_expression
, test_and_plot
)
## #### ####
## Distribution tests
## #### ####
## # A tibble: 4 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 8 BGLU2 0.734 0.672 0.579
## 2 8 HSP70 0.734 0.672 0.579
## 3 8 PTI5 0.734 0.672 0.579
## 4 8 13-LOX 0.937 0.921 0.788
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 0.752 0.419
## 2 HSP70 1 6 1.64 0.248
## 3 PTI5 1 6 0.310 0.598
## 4 13-LOX 1 6 1.36 0.288
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 0.752 0.419
## 2 HSP70 1 6 1.64 0.248
## 3 PTI5 1 6 0.310 0.598
## 4 13-LOX 1 6 1.36 0.288
## #### ####
## Fligner
## #### ####
## # A tibble: 4 × 5
## # Groups: transcript [4]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 BGLU2 1.75 0.186 1 Fligner-Killeen test of homogeneity of…
## 2 HSP70 0.715 0.398 1 Fligner-Killeen test of homogeneity of…
## 3 PTI5 0.306 0.580 1 Fligner-Killeen test of homogeneity of…
## 4 13-LOX 0.619 0.432 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.102 BGLU2 4 4 small
## 2 measurement noninoculated inoculated 0.408 HSP70 4 4 moderate
## 3 measurement noninoculated inoculated 0.612 PTI5 4 4 large
## 4 measurement noninoculated inoculated 0.204 13-LOX 4 4 small
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated -0.0596 BGLU2 4 4 negligible
## 2 measurement noninoculated inoculated -0.800 HSP70 4 4 moderate
## 3 measurement noninoculated inoculated -1.16 PTI5 4 4 large
## 4 measurement noninoculated inoculated -0.497 13-LOX 4 4 small
## [1] transcript group1 group2 p p.adj
## [6] p.adj.signif
## <0 rows> (or 0-length row.names)
results_1h = run_analysis_for_time(data = PS218
, time_point = "1 h"
, myorder
, pal
, groupvars
, measurevar
, y_scales1
, y_scales2
, my_dir
, plot_gene_expression
, test_and_plot
)
## #### ####
## Distribution tests
## #### ####
## # A tibble: 4 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 7 BGLU2 0.363 NA 0.718
## 2 7 HSP70 0.0973 NA 0.188
## 3 7 PTI5 0.0782 NA 0.188
## 4 7 13-LOX 0.363 NA 0.718
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 5 93.6 0.000200
## 2 HSP70 1 5 5.10 0.0736
## 3 PTI5 1 5 1.07 0.348
## 4 13-LOX 1 5 1.85 0.232
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 5 93.6 0.000200
## 2 HSP70 1 5 5.10 0.0736
## 3 PTI5 1 5 1.07 0.348
## 4 13-LOX 1 5 1.85 0.232
## #### ####
## Fligner
## #### ####
## # A tibble: 4 × 5
## # Groups: transcript [4]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 BGLU2 4.01 0.0453 1 Fligner-Killeen test of homogeneity of…
## 2 HSP70 4.05 0.0441 1 Fligner-Killeen test of homogeneity of…
## 3 PTI5 0.662 0.416 1 Fligner-Killeen test of homogeneity of…
## 4 13-LOX 1.89 0.170 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0 BGLU2 3 4 small
## 2 measurement noninoculated inoculated 0.401 HSP70 3 4 moderate
## 3 measurement noninoculated inoculated 0.668 PTI5 3 4 large
## 4 measurement noninoculated inoculated 0.802 13-LOX 3 4 large
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated -0.393 BGLU2 3 4 small
## 2 measurement noninoculated inoculated -0.912 HSP70 3 4 large
## 3 measurement noninoculated inoculated -1.35 PTI5 3 4 large
## 4 measurement noninoculated inoculated -3.09 13-LOX 3 4 large
## [1] transcript group1 group2 p p.adj
## [6] p.adj.signif
## <0 rows> (or 0-length row.names)
results_2h = run_analysis_for_time(data = PS218
, time_point = "2 h"
, myorder
, pal
, groupvars
, measurevar
, y_scales1
, y_scales2
, my_dir
, plot_gene_expression
, test_and_plot
)
## #### ####
## Distribution tests
## #### ####
## # A tibble: 4 × 7
## n transcript Shapiro_Wilk_BH Anderson_Darling_BH Lilliefors_KS_BH
## <int> <fct> <dbl> <dbl> <dbl>
## 1 8 BGLU2 0.215 0.290 0.755
## 2 8 HSP70 0.0287 0.0275 0.0291
## 3 8 PTI5 0.0959 0.127 0.400
## 4 8 13-LOX 0.0287 0.0275 0.0914
## # ℹ 2 more variables: Jarque_Bera_BH <dbl>, DAgostino_Skewness_BH <dbl>
## #### ####
## Quantile-Quantile plots
## #### ####
## #### ####
## Test for homogeneity of variance across groups
## #### ####
## #### ####
## Levene
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 0.151 0.711
## 2 HSP70 1 6 0.926 0.373
## 3 PTI5 1 6 2.14 0.194
## 4 13-LOX 1 6 0.634 0.456
## #### ####
## Brown-Forsythe
## #### ####
## # A tibble: 4 × 5
## transcript df1 df2 statistic p
## <fct> <int> <int> <dbl> <dbl>
## 1 BGLU2 1 6 0.151 0.711
## 2 HSP70 1 6 0.926 0.373
## 3 PTI5 1 6 2.14 0.194
## 4 13-LOX 1 6 0.634 0.456
## #### ####
## Fligner
## #### ####
## # A tibble: 4 × 5
## # Groups: transcript [4]
## transcript statistic p.value parameter method
## <fct> <dbl> <dbl> <dbl> <chr>
## 1 BGLU2 0.0945 0.759 1 Fligner-Killeen test of homogeneity of…
## 2 HSP70 0.120 0.729 1 Fligner-Killeen test of homogeneity of…
## 3 PTI5 1.45 0.228 1 Fligner-Killeen test of homogeneity of…
## 4 13-LOX 1.46 0.226 1 Fligner-Killeen test of homogeneity of…
## #### ####
## Wilcoxon effect size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.102 BGLU2 4 4 small
## 2 measurement noninoculated inoculated 0 HSP70 4 4 small
## 3 measurement noninoculated inoculated 0.816 PTI5 4 4 large
## 4 measurement noninoculated inoculated 0.816 13-LOX 4 4 large
## #### ####
## Cohen's d Measure of Effect Size
## #### ####
## # A tibble: 4 × 8
## .y. group1 group2 effsize transcript n1 n2 magnitude
## * <chr> <chr> <chr> <dbl> <fct> <int> <int> <ord>
## 1 measurement noninoculated inoculated 0.413 BGLU2 4 4 small
## 2 measurement noninoculated inoculated -0.442 HSP70 4 4 small
## 3 measurement noninoculated inoculated -2.37 PTI5 4 4 large
## 4 measurement noninoculated inoculated -12.6 13-LOX 4 4 large
## transcript group1 group2 p p.adj p.adj.signif
## 3 PTI5 noninoculated inoculated 0.003937008 0.01574803 *
## 4 13-LOX noninoculated inoculated 0.023622047 0.04724409 *
7 sessionInfo
## R version 4.4.1 (2024-06-14 ucrt)
## Platform: x86_64-w64-mingw32/x64
## Running under: Windows 11 x64 (build 26100)
##
## Matrix products: default
##
##
## locale:
## [1] LC_COLLATE=English_United Kingdom.utf8
## [2] LC_CTYPE=English_United Kingdom.utf8
## [3] LC_MONETARY=English_United Kingdom.utf8
## [4] LC_NUMERIC=C
## [5] LC_TIME=English_United Kingdom.utf8
##
## time zone: Europe/Ljubljana
## tzcode source: internal
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] ggplot2_4.0.0 magrittr_2.0.3
##
## loaded via a namespace (and not attached):
## [1] DBI_1.2.3 Rdpack_2.6.4 sandwich_3.1-1
## [4] MKinfer_1.2 rlang_1.1.5 multcomp_1.4-28
## [7] miceadds_3.18-36 tseries_0.10-58 matrixStats_1.5.0
## [10] compiler_4.4.1 vctrs_0.6.5 quadprog_1.5-8
## [13] pkgconfig_2.0.3 shape_1.4.6.1 crayon_1.5.3
## [16] fastmap_1.2.0 backports_1.5.0 labeling_0.4.3
## [19] utf8_1.2.6 rmarkdown_2.30 exactRankTests_0.8-35
## [22] nloptr_2.2.1 purrr_1.0.4 xfun_0.53
## [25] glmnet_4.1-10 modeltools_0.2-24 jomo_2.7-6
## [28] cachem_1.1.0 jsonlite_2.0.0 gmp_0.7-5
## [31] pan_1.9 broom_1.0.10 parallel_4.4.1
## [34] R6_2.6.1 coin_1.4-3 bslib_0.9.0
## [37] stringi_1.8.7 RColorBrewer_1.1-3 rpart_4.1.24
## [40] car_3.1-3 boot_1.3-31 jquerylib_0.1.4
## [43] Rcpp_1.0.14 iterators_1.0.14 knitr_1.50
## [46] zoo_1.8-14 nnet_7.3-20 Matrix_1.7-1
## [49] splines_4.4.1 tidyselect_1.2.1 rstudioapi_0.17.1
## [52] dichromat_2.0-0.1 abind_1.4-8 yaml_2.3.10
## [55] codetools_0.2-20 curl_6.2.2 arrangements_1.1.9
## [58] lattice_0.22-6 tibble_3.3.0 quantmod_0.4.28
## [61] withr_3.0.2 S7_0.2.0 evaluate_1.0.5
## [64] moments_0.14.1 survival_3.8-3 zip_2.3.3
## [67] xts_0.14.1 pillar_1.11.1 ggpubr_0.6.1
## [70] carData_3.0-5 mice_3.18.0 nortest_1.0-4
## [73] foreach_1.5.2 stats4_4.4.1 reformulas_0.4.1
## [76] generics_0.1.4 TTR_0.24.4 scales_1.4.0
## [79] minqa_1.2.8 glue_1.8.0 tools_4.4.1
## [82] data.table_1.17.8 MKdescr_0.9 lme4_1.1-37
## [85] openxlsx_4.2.8 ggsignif_0.6.4 mvtnorm_1.3-3
## [88] cowplot_1.2.0 grid_4.4.1 mitools_2.4
## [91] tidyr_1.3.1 rbibutils_2.3 libcoin_1.0-10
## [94] nlme_3.1-166 Formula_1.2-5 cli_3.6.3
## [97] dplyr_1.1.4 ggh4x_0.3.1 gtable_0.3.6
## [100] rstatix_0.7.2 sass_0.4.10 digest_0.6.37
## [103] TH.data_1.1-4 farver_2.1.2 htmltools_0.5.8.1
## [106] lifecycle_1.0.4 mitml_0.4-5 MASS_7.3-64